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ABSTRACT 

We  investigate  the  use  of  optical  photometric  variability  to  select  and  identify  blazars  in  large-scale  time-domain 
surveys,  in  part  to  aid  in  the  identification  of  blazar  counterparts  to  the  ~30%  of  y- ray  sources  in  the  Fermi  2FGL 
catalog  still  lacking  reliable  associations.  Using  data  from  the  optical  LINEAR  asteroid  survey,  we  characterize  the 
optical  variability  of  blazars  by  fitting  a  damped  random  walk  model  to  individual  light  curves  with  two  main  model 
parameters,  the  characteristic  timescales  of  variability  r,  and  driving  amplitudes  on  short  timescales  a .  Imposing 
cuts  on  minimum  t  and  a  allows  for  blazar  selection  with  high  efficiency  E  and  completeness  C.  To  test  the 
efficacy  of  this  approach,  we  apply  this  method  to  optically  variable  LINEAR  objects  that  fall  within  the  several- 
arcminute  error  ellipses  of  y- ray  sources  in  the  Fermi  2FGL  catalog.  Despite  the  extreme  stellar  contamination  at 
the  shallow  depth  of  the  LINEAR  survey,  we  are  able  to  recover  previously  associated  optical  counterparts  to  Fermi 
active  galactic  nuclei  with  E  ^  88%  and  C  =  88%  in  Fermi  95%  confidence  error  ellipses  having  semimajor  axis 
r  <  8'.  We  find  that  the  suggested  radio  counterpart  to  Fermi  source  2FGL  J1649. 6+5238  has  optical  variability 
consistent  with  other  y- ray  blazars  and  is  likely  to  be  the  y- ray  source.  Our  results  suggest  that  the  variability  of  the 
non-thermal  jet  emission  in  blazars  is  stochastic  in  nature,  with  unique  variability  properties  due  to  the  effects  of 
relativistic  beaming.  After  correcting  for  beaming,  we  estimate  that  the  characteristic  timescale  of  blazar  variability 
is  ~3  years  in  the  rest  frame  of  the  jet,  in  contrast  with  the  ~320  day  disk  flux  timescale  observed  in  quasars.  The 
variability-based  selection  method  presented  will  be  useful  for  blazar  identification  in  time-domain  optical  surveys 
and  is  also  a  probe  of  jet  physics. 
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1.  INTRODUCTION 

Blazars  are  a  relatively  rare  sub-class  of  active  galactic 
nuclei  (AGNs)  in  which  a  jet  is  aligned  along  the  observer’s 
line  of  sight,  leading  to  the  effects  of  relativistic  beaming  and 
unusual  associated  emission  (Blandford  &  Rees  1978).  Blazars 
are  among  the  most  variable  extragalactic  objects  detected  in 
time-domain  optical  surveys  and  have  strong  emission  from 
radio  to  TeV  energies  (Ulrich  et  al.  1997).  The  central  engine  is 
believed  to  be  accretion  onto  a  supermassive  black  hole,  driving 
relativistic  outflows  in  a  collimated  jet  with  typical  Lorentz 
factors  on  the  order  of  T  ~  10  that  is  pointed  to  within  angle 
T_1  of  the  observer.  The  term  “blazars”  usually  encompasses 
both  BL  Lac  objects  and  flat-spectrum  radio  quasars  (FSRQs), 
which  are  believed  to  be  jet- aligned  Fanaroff-Riley  type  I  and 
II  AGNs,  respectively  (Urry  &  Padovani  1995).  In  this  paper, 
we  adopt  this  definition. 

The  canonical  broadband  spectral  energy  distribution  (SED) 
of  blazars  typically  includes  several  main  components:  (1)  a 
synchrotron  peak  likely  due  to  tangled  magnetic  fields  in  the  jet 
that  may  extend  from  the  radio  to  the  soft-X-ray  regime,  (2)  an 
inverse-Compton  peak  in  the  X-ray  to  GeV  regime  likely  due  to 
scattering  of  synchrotron  or  external  photons  off  of  relativistic 
electrons  in  the  jet,  (3)  a  possible  inverse-Compton  component 


in  the  soft  X-rays  due  to  a  hot  corona,  and  (4)  occasional  hints  of 
an  underlying  accretion  disk  continuum  or  host  galaxy  emission. 
The  optical  and  y- ray  observations  we  use  in  this  paper  are 
expected  to  be  dominated  by  the  synchrotron  and  jet  inverse- 
Compton  emission,  respectively. 

The  strong  high-energy  inverse-Compton  peak  in  the  SEDs 
of  blazars  causes  them  to  account  for  the  vast  majority  of 
bright  extragalactic  y -ray-emitting  sources.  The  Fermi  Space 
Telescope  has  surveyed  the  y- ray  sky  since  launch  in  2008,  and 
its  Large  Area  Telescope  (LAT;  Atwood  et  al.  2009)  instrument 
provides  by  far  the  deepest  survey  to  date  in  the  100  MeV  to 
100  GeV  regime.  The  current  2  year  LAT  source  catalog  (2FGL; 
Nolan  et  al.  2012)  includes  1873  total  sources,  ~44%  of  which 
are  reliably  associated  with  AGNs,  and  an  additional  M4%  are 
candidate  AGN  associations.  Of  the  reliably  associated  AGNs, 
the  overwhelming  majority  are  blazars.  Approximately  half  of 
these  Fermi  blazars  are  BL  Lac  objects,  and  half  are  FSRQs 
(Ackermann  et  al.  201 1). 

The  2FGL  catalog  is  a  significant  improvement  over  the  1 1 
month  Fermi  LAT  source  catalog  (Abdo  et  al.  2010)  in  num¬ 
ber  of  sources,  source  detection  methods,  and  source  associ¬ 
ations  to  known  objects.  However,  ~32%  of  sources  in  the 
2FGL  catalog  still  lack  reliable  associations  with  any  known  ob¬ 
ject,  many  of  which  may  be  blazars,  especially  at  high  Galactic 
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latitudes  (Ackermann  et  al.  2012).  Although  the  angular  resolu¬ 
tion  of  the  LAT  is  a  dramatic  improvement  over  previous  all- sky 
]/-ray  survey  instruments  such  as  EGRET,  typical  95%  confi¬ 
dence  error  ellipses  of  sources  in  the  2FGL  catalog  are  still  on 
the  order  of  several  arcminutes  in  size  and  can  contain  numerous 
candidate  counterparts  at  other  wavelengths,  making  association 
and  identification  of  y-ray  sources  difficult. 

Recently,  huge  efforts  have  gone  into  identifying  Fermi 
y- ray  sources  by  positional  coincidences  of  candidate  counter¬ 
parts  from  multi- wavelength  surveys  (e.g.,  Stephen  et  al.  2010; 
Maeda  et  al.  2011),  observations  of  contemporaneous  variabil¬ 
ity  at  different  wavelengths  (e.g.,  Kara  et  al.  2012),  and  statis¬ 
tical  methods  based  on  observed  y-ray  source  properties  (e.g., 
Ackermann  et  al.  2012).  Notably,  Massaro  et  al.  (2011)  ap¬ 
plied  a  variant  of  the  well-known  method  of  selecting  AGNs  by 
their  unique  colors  in  the  mid-IR  (e.g.,  Lacy  et  al.  2004;  Stem 
et  al.  2005)  to  data  from  the  WISE  survey  (Wright  et  al.  2010) 
and  used  their  selection  method  to  separate  blazars  from  stars. 
Massaro  et  al.  (2012)  then  further  utilized  their  technique  to  find 
WISE- detected  y-ray  blazar  candidates  in  the  2FGL,  achieving 
excellent  completeness,  although  the  efficiency  of  their  method 
is  unclear  (see  discussion  in  Section  5.1). 

Aside  from  their  distinct  mid-IR  colors,  blazars  are  also 
unique  in  their  strong  variability  at  nearly  all  wavelength 
regimes.  This  will  make  blazars  stand  out  in  the  flood  of 
time-domain  data  from  current  and  future  large-scale  optical 
imaging  surveys  such  as  Pan-STARRS  (Kaiser  et  al.  2010), 
the  Catalina  Real-Time  Transient  Survey  (Drake  et  al.  2009), 
the  Palomar  Transient  Factory  (PTF;  Rau  et  al.  2009),  and  the 
Large  Synoptic  Survey  Telescope  (LSST;  Ivezic  et  al.  2008).  It 
is  thus  an  auspicious  time  to  explore  the  possibility  of  selecting 
and  identifying  y  -ray-emitting  blazars  en  masse  in  Fermi  error 
ellipses  by  their  optical  photometric  variability. 

Blazars  have  observed  variability  on  timescales  from  hours  to 
years,  and  their  distinctiveness  in  time-domain  surveys  results 
from  two  effects.  First,  the  relativistic  jet  strongly  beams  non- 
thermal  emission  along  the  line  of  sight  toward  the  observer, 
boosting  the  luminosity  by  several  orders  of  magnitude.  Second, 
relativistic  Doppler-boosting  shortens  the  observed  timescale  of 
variability  in  comparison  to  that  in  the  rest  frame  of  the  outflow. 
However,  the  physical  mechanisms  causing  the  variability  are 
far  less  clear.  Since  blazar  emission  is  largely  dominated  by  non- 
thermal  emission  from  the  jet,  it  is  possible  that  internal  shocks 
from  overtaking  collisions  of  fluid  shells  with  different  velocities 
along  the  jet  produce  time-variable  synchrotron  emission  (e.g., 
Bottcher  &  Dermer  2010),  although  other  emission  components 
such  as  the  accretion  disk  may  also  be  time  variable. 

Other  classes  of  AGNs,  including  the  far  more  numerous  “nor¬ 
mal”  quasars  (i.e.,  non-FSRQ,  type  1  AGNs),  are  also  highly  op¬ 
tically  variable,  although  the  predominance  of  particular  physi¬ 
cal  variability  mechanisms  may  differ  among  AGN  subclasses. 
Sesar  et  al.  (2007)  have  shown  that  ~90%  of  quasars  in  Sloan 
Digital  Sky  Survey  (SDSS;  York  et  al.  2000)  Stripe  82  are 
variable  at  the  >0.03  intrinsic  rms  variability  level  (defined  in 
Equation  (1)),  and  previous  studies  using  large  samples  of  indi¬ 
vidual  quasar  light  curves  have  shown  that  quasar  optical  vari¬ 
ability  is  well  described  by  a  damped  random  walk  (DRW) 
model  (Kelly  et  al.  2009;  MacLeod  et  al.  2010;  Kozlowski  et  al. 
2010;  Butler  &  Bloom  2011;  Zu  et  al.  2012).  MacLeod  et  al. 
(2011,  hereafter  MA11)  further  utilized  the  DRW  model  as  a 
tool  to  separate  stars  and  quasars  by  optical  photometric  vari¬ 
ability  in  SDSS  Stripe  82  and  showed  that  quasar  selection  using 
this  method  can  achieve  efficiency  E  ^  80%  and  completeness 


C  =  90%  by  variability  alone.  This  is  particularly  useful  for 
selection  of  quasars  in  certain  redshift  regimes,  where  color  se¬ 
lection  may  fail  due  to  overwhelming  contamination  from  fore¬ 
ground  stars.  Furthermore,  studies  of  the  long-term  variability 
characteristics  of  AGNs  may  also  provide  a  unique  perspective 
on  accretion  disk  and  jet  physics. 

The  optical  variability  of  small  numbers  of  individual  blazars 
has  been  studied  extensively  in  the  literature  through  inten¬ 
sive  long-term  monitoring  campaigns  (e.g.,  Webb  et  al.  1988; 
Carini  et  al.  1992).  However,  no  consistent  picture  of  variabil¬ 
ity  has  emerged,  at  least  partially  due  to  the  heterogenous  na¬ 
ture  of  these  observations  and  the  small  sample  sizes.  Blazar 
light  curves  from  a  large-scale  flux-limited  time-domain  survey 
should  be  more  homogenous,  and  the  depth  of  current  surveys 
can  yield  orders  of  magnitude  more  blazar  light  curves  than 
possible  through  targeted  monitoring  campaigns.  However,  the 
rarity  of  blazars  places  significant  constraints  on  the  specifi¬ 
cations  of  any  time-domain  survey  in  which  a  large  number 
of  blazar  light  curves  can  be  obtained.  Despite  the  success  of 
the  SDSS  Stripe  82  for  time-domain  astronomy,  the  depth  of  the 
survey  (~22.5  mag  in  r  band)  does  not  compensate  for  the  small 
sky  coverage  (~290  deg2)  when  it  comes  to  studying  blazars. 
We  instead  employ  the  recalibrated  LINEAR  survey,  which  es¬ 
sentially  covers  the  ~  10,000  deg2  SDSS  photometric  footprint 
down  to  ~17  mag  in  r  band  (Sesar  et  al.  2011,  hereafter  SE1 1). 
Bauer  et  al.  (2009a)  provided  a  previous  study  of  blazar  optical 
variability  in  the  Palomar-Quest  Survey,  and  Bauer  et  al.  (2009b) 
conducted  a  subsequent  variability-based  search  for  new  blazars 
within  that  survey.  However,  the  Palomar-Quest  survey  provided 
a  median  of  ~5  epochs  of  observation  for  each  object  per  fil¬ 
ter  and  could  only  be  used  to  study  the  ensemble  variability  of 
blazars  rather  than  the  individual  sources.  By  contrast,  the  ~200 
epochs  of  observation  per  object  in  the  LINEAR  survey  allow 
us  to  directly  model  the  light  curves  of  individual  blazars  in  this 
study. 

The  structure  of  this  paper  is  as  follows:  In  Section  2, 
we  describe  the  LINEAR  survey  and  the  construction  of 
the  Variable  LINEAR  catalog,  from  which  our  blazar  and 
comparison  normal  quasar  samples  are  drawn.  In  Section  3, 
we  describe  the  DRW  model  of  variability  and  our  results  of 
modeling  the  LINEAR  light  curves  of  blazars,  normal  quasars, 
and  stars.  In  Section  4,  we  test  the  efficiency  and  completeness  of 
blazar  selection  using  a  DRW  model  within  Fermi  error  ellipses. 
In  Section  5,  we  discuss  this  method  in  the  context  of  its  future 
applications,  as  well  as  possible  implications  of  our  work  for 
AGN  physics.  We  summarize  and  conclude  in  Section  6. 

2.  THE  LINEAR  SURVEY 
2.7.  Photometric  Data 

All  photometric  data  used  in  this  paper  are  from  the  archives 
of  the  MIT  Lincoln  Laboratory  (MITLL)  Lincoln  Near-Earth 
Asteroid  Research  (LINEAR)  survey,  spanning  the  period  from 
2002  December  through  2008  March.  A  review  of  the  original 
LINEAR  near-Earth  asteroid  survey  program  is  presented  in 
Stokes  et  al.  (2000),  and  the  subsequent  photometric  recalibra¬ 
tion  of  the  archived  data  using  SDSS  to  construct  the  LINEAR 
photometric  database  is  discussed  in  SE1 1.  We  summarize  only 
the  most  salient  points  here.  The  LINEAR  survey  program  used 
two  1.01  m  diameter  telescopes  at  the  Experimental  Test  Site 
within  the  US  Army  White  Sands  Missile  Range  in  New  Mex¬ 
ico,  each  equipped  with  a  5  megapixel  (2560  x  1960)  back- 
illuminated  CCD  developed  at  MITLL  and  described  in  Burke 
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et  al.  (1998).  The  recalibrated  survey  covers  ~  10,000  deg2  of 
sky,  overlapping  the  SDSS  photometric  footprint,  and  contains 
over  5  billion  photometric  measurements  of  ~25  million  objects 
down  to  a  5a  depth  of  r  <  18  mag.  Due  to  the  astrometric  goals 
of  the  original  survey  and  to  increase  signal-to-noise  ratio,  only 
a  single  broad  filter  was  used.  The  cadence  range  spans  minutes 
to  years,  with  a  peak  at  the  main  15  minute  cadence,  a  gap  at 
8  hr,  and  a  secondary  peak  around  1 1  days.  The  median  number 
of  good  observations  per  object  in  the  full  LINEAR  catalog  is 
~460  within  ±10°  of  the  ecliptic  plane  and  ~200  elsewhere. 

SE11  extracted  sources  from  LINEAR  imaging  using  fixed- 
aperture  photometry  and  recalibrated  the  astrometry  to  the 
USNO-B  catalog  (Barron  et  al.  2008).  To  recalibrate  the  pho¬ 
tometry  using  SDSS,  SE11  first  matched  LINEAR  sources  to 
SDSS  DR7  non-saturated,  primary  objects.  Using  gri  photome¬ 
try  of  DR7  objects  matched  to  LINEAR,  SE1 1  then  modeled  and 
calculated  synthetic  LINEAR  magnitudes  of  matched  SDSS  ob¬ 
jects  rasDSS,  effectively  turning  the  SDSS  imaging  catalog  into 
a  catalog  of  LINEAR  photometric  standard  stars.  After  a  super 
flat-field  correction  on  each  field  using  these  calibration  stars, 
a  catalog  of  5  billion  individual  point  sources  with  recalibrated 
LINEAR  magnitudes  m linear  from  all  epochs  of  observations 
was  positionally  clustered  into  25  million  objects,  each  with  at 
least  15  epochs  of  observation.  This  constitutes  the  full  recali¬ 
brated  LINEAR  catalog.  Checks  on  the  recalibrated  photometry 
by  SE11  show  that  the  median  w linear  —  ^sdss  residual  per 
field  has  a  distribution  of  about  0.01  mag  wide.  This  should  not 
be  confused  with  the  single-epoch  photometric  uncertainties  for 
individual  objects  in  LINEAR,  which  are  generally  ->0.04  mag 
(see  Figure  2).  Fields  with  w linear  —  msDss  >0.1  (usually  due 
to  variable  cloud  coverage)  are  removed. 

2.2.  The  Bright  LINEAR  Catalog 

From  the  full  LINEAR  catalog,  we  select  a  bright  subset 
catalog  suitable  for  variability  science  by  imposing  a  cut  on  the 
minimum  number  of  good  observations  in  each  light  curve  of 
^30  and  a  LINEAR  recalibrated  magnitude  cut  of  14  <  iulinear 
<17  (where  photometric  errors  are  <0.1 1  mag).  At  magnitudes 
>  17  and  <14,  single-epoch  photometric  errors  rise  rapidly  due 
to  photon  noise  and  saturation,  respectively.  These  cuts  provide 
us  with  ~4.5  million  objects  in  the  Bright  LINEAR  catalog. 

To  create  a  sample  of  known  blazars  and  quasars  in  the  Bright 
LINEAR  catalog,  we  positionally  match  the  LINEAR  sources  to 
catalogs  of  known  BL  Lac  objects  and  FSRQs.  To  find  the  known 
BL  Lac  objects,  we  match  to  1371  BL  Lac  objects  in  the  Veron- 
Cetty  &  Veron  (2010)  catalog,  501  BL  Lac  objects  in  the  Plotkin 
et  al.  (2008)  catalog,  and  637  radio-loud  BL  Lac  objects  in  the 
Plotkin  et  al.  (2010)  catalog  (none  of  the  86  radio-quiet  objects 
in  Plotkin  et  al.  (2010)  matched  to  a  LINEAR  object),  all  using 
a  3"0  matching  radius,  resulting  in  101  matched  distinct  BL  Lac 
objects.  Known  quasars  in  the  Bright  LINEAR  catalog  were 
identified  by  matching  to  the  105,782  quasars  in  the  Schneider 
et  al.  (2010)  SDSS  Data  Release  7  catalog,  also  with  a  3"0 
matching  radius,  resulting  in  1020  matched  distinct  quasars. 
We  also  separate  out  the  FSRQs  by  positionally  matching  these 
1020  quasars  to  the  CRATES  radio  survey  (Healey  et  al.  2007) 
of  11,131  flat-spectrum  (a  >  —0.5,  where  Sv  oc  v01)  radio 
sources  using  a  3'.'0  matching  radius,  resulting  in  42  distinct 
FSRQs.  The  978  non-flat  radio  spectrum  quasars  remaining  are 
classified  here  as  normal  quasars  (i.e.,  unlikely  to  be  strongly 
dominated  in  the  optical  by  jet  emission). 

We  note  that  there  is  ultimately  an  overlap  of  three  objects 
between  the  BL  Lac  object  and  FSRQ  samples  due  to  double 


Figure  1.  Normalized  histograms  of  the  number  of  light-curve  epochs  available 
for  each  source  in  the  Bright  LINEAR  catalog.  The  distributions  are  shown  for 
140  blazars  (solid),  978  normal  quasar  (dashed),  and  all  other  objects  (dotted). 

counting  objects  that  appear  in  the  BL  Lac  object  catalogs, 
the  Schneider  et  al.  (2010)  quasar  catalog,  and  the  Healey 
et  al.  (2007)  catalog,  which  may  result  from  the  intermittent 
appearance  of  broad  emission  lines  in  their  optical  spectra.  Since 
our  definition  of  “blazars”  here  includes  both  BL  Lac  objects 
and  FSRQs,  we  count  these  three  objects  as  blazars  (but  not 
normal  quasars  due  to  their  flat  radio  spectra).  In  summary,  we 
identified  a  total  of  140  blazars  and  978  normal  quasars  in  the 
Bright  LINEAR  catalog. 

The  potential  of  the  LINEAR  survey  for  blazar  time-domain 
studies  can  be  seen  in  Figure  1,  which  shows  normalized 
histograms  of  the  number  of  LINEAR  epochs  of  observation 
for  each  of  the  140  objects  in  the  blazar  sample,  the  987  objects 
in  the  normal  quasar  sample,  and  all  other  objects  (mainly 
stars)  in  the  Bright  LINEAR  catalog.  The  histograms  are  all 
generally  consistent  with  each  other,  peak  at  ~200  epochs,  and 
lack  objects  with  <30  good  observations  due  to  the  previously 
imposed  cut. 

Before  attempting  to  model  the  light  curve  of  every  object, 
we  would  first  like  to  compare  the  general  level  of  variability  of 
blazars,  normal  quasars,  and  other  objects  in  the  Bright  LINEAR 
catalog.  We  follow  Sesar  et  al.  (2007)  and  estimate  the  intrinsic 
rms  variability  a  for  each  light  curve,  defined  as 

=  P?  —  iA(^linear)2]1/2  (1) 

for  all  objects  with  X  >  \l/(m linear)  and  a  =  0  otherwise, 
where  X  is  the  rms  scatter  of  m linear  for  all  observations  of 
each  object,  and  linear)  is  the  median  photometric  error  of 
LINEAR  objects  as  a  function  of  magnitude  found  by  SE1 1  (see 
their  Equation  (7)).  Since  the  photometric  errors  in  LINEAR  are 
relatively  large,  X  mainly  reflects  photometric  noise  rather  than 
variability  for  the  vast  majority  of  objects.  Equation  (1)  attempts 
to  remove  the  effects  of  the  median  photometric  error,  and  thus 
a  is  a  better  measure  of  intrinsic  variability.  Figure  2  shows  the 
median  a  and  \l/(m linear)  of  all  objects  in  the  Bright  LINEAR 
catalog  as  a  function  of  wlinear  magnitude. 

We  compare  the  variability  of  blazars  to  normal  quasars  and 
other  LINEAR  objects,  as  it  is  expected  that  blazars  should 
be  systematically  more  variable.  Figure  3  shows  the  integrated 
distribution  of  blazars,  normal  quasars,  and  other  LINEAR 
bright  objects  in  the  Bright  LINEAR  catalog  as  a  function  of  the 
estimated  intrinsic  rms  variability  a.  We  estimate  that  ~36% 
of  blazars,  as  compared  to  only  ~11%  of  quasars,  are  variable 
with  a  >  0.15  mag  in  this  passband.  At  this  level  of  variability, 
the  estimate  of  a  is  not  significantly  affected  by  uncertainties  in 
VK^linear)-  Only  ~4%  of  the  other  bright  LINEAR  objects 
(non-blazars  and  non-quasars)  are  variable  above  this  a  > 
0.15  mag  level,  and  they  are  likely  a  combination  of  variable 
stars  (e.g.,  RR  Lyrae),  eclipsing  binaries,  underestimates  of  the 
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Figure  2.  Median  intrinsic  rms  variability  a  (solid)  and  photometric  uncertainty 
VK^linear)  (dash-dotted;  Equation  (7)  from  SE11)  of  LINEAR  objects  in  the 
Bright  LINEAR  catalog,  as  a  function  of  recalibrated  LINEAR  magnitude. 

photometric  errors  x/r,  source  blending,  or  other  photometric 
issues. 

2.3.  The  Variable  LINEAR  Catalog 

In  order  to  characterize  the  photometric  variability,  we  need 
a  clean  sample  of  light  curves  dominated  by  intrinsic  variability 
rather  than  photometric  errors.  Thus,  we  follow  the  method  of 
Sesar  et  al.  (2007)  to  create  a  clean  Variable  LINEAR  catalog 
as  a  subsample  of  the  4.5  million  objects  in  the  Bright  LINEAR 
catalog  discussed  in  the  previous  section. 

As  a  first  variability  criterion,  we  require  a  >  0.1,  roughly 
equivalent  to  selecting  sources  with  a  more  than  two  times 
the  measurement  noise  at  14  <  m  linear  <16,  and  roughly 
corresponding  to  the  photometric  error  at  m  linear  ~  17. 
Approximately  8%  of  objects  in  the  Bright  LINEAR  catalog 
pass  this  initial  selection  cut,  including  some  non-variable 
objects  at  the  faint  end  that  have  large  a  due  to  large  photometric 
errors.  To  reduce  these  spurious  contaminants,  we  assume  that 
the  photometric  error  distribution  follows  a  Gaussian  and  place 
a  cut  on  the  x  2  per  degree  of  freedom  calculated  with  respect  to  a 
weighted  mean  magnitude  and  errors  from  the  photometry.  A  cut 
of  x2/dof  >  3  leaves  188,745  objects  composing  the  Variable 
LINEAR  catalog  we  use  hereafter,  each  typically  having  ~200 
epochs  per  light  curve. 

Matching  this  variable  LINEAR  sample  to  the  various  cata¬ 
logs  described  in  Section  2.2  yields  60  blazars  and  155  normal 
quasars,  including  one  of  the  three  overlapping  objects  discussed 
previously  in  Section  2.2  that  we  include  in  both  samples.  Al¬ 
though  the  variable  LINEAR  sample  cuts  out  the  majority  of 
the  known  140  blazars  and  978  normal  quasars  in  the  Bright 
LINEAR  catalog,  this  is  not  because  the  majority  of  blazars  and 
normal  quasars  are  non-variable,  but  rather  because  the  large 
photometric  errors  in  LINEAR  overwhelm  the  intrinsic  vari¬ 
ability  of  fainter  LINEAR  objects.  This  effect  can  be  seen  in  the 
much  larger  fraction  of  blazars  (60  out  of  140)  that  survive  these 
variability  cuts,  in  comparison  with  normal  quasars  (155  out  of 
978).  Since  blazars  tend  to  be  systematically  more  variable  than 
normal  quasars  (as  we  have  shown  in  Section  2.2),  their  intrin¬ 
sic  variability  will  dominate  over  photometric  errors  for  a  much 
larger  fraction  of  objects. 

3.  THE  DAMPED  RANDOM  WALK  MODEL 

MAI  1  designed  and  tested  optical  variability-based  selection 
of  quasars  in  SDSS  Stripe  82  by  modeling  the  individual  light 
curves  of  known  quasars  as  a  DRW.  We  perform  a  similar 
analysis  of  the  LINEAR  light  curves  of  known  blazars  and 
normal  quasars  and  show  that  blazar  light  curves  lie  in  a  distinct 
region  of  DRW  variability  parameter  space.  The  DRW  model 


a  [mag] 


Figure  3.  Fraction  of  objects  with  an  intrinsic  rms  variability  larger  than  cr,  as 
a  function  of  cr,  for  blazars  (solid),  normal  quasars  (dashed),  and  other  objects 
(dash-dotted). 

statistically  parameterizes  variability  using  three  parameters:  a 
mean  light-curve  magnitude,  a  characteristic  damping  timescale 
of  variability  r,  and  a  driving  amplitude  of  short-term  stochastic 
variability  o  (see  Kelly  et  al.  2009;  Kozlowski  et  al.  2010; 
MacLeod  et  al.  2010;  Zu  et  al.  2012,  for  further  discussion). 
The  fitted  parameters  can  be  expressed  in  terms  of  the  structure 
function  (SF),  defined  as  the  rms  magnitude  difference  between 
all  pairs  of  observations  in  each  individual  light  curve  as  a 
function  of  the  time  lag.  While  this  model  is  phenomenological, 
it  is  likely  that  the  fitted  parameters  reflect  the  physical  processes 
that  cause  the  variability. 

Following  Kozlowski  et  al.  (2010)  and  MA11,  we  model 
each  individual  light  curve  in  our  sample  as  a  stochastic  process 
described  by  the  exponential  covariance  matrix 

Sij  =  o-2exp(— I ti  -  tj\/x)  (2) 

for  each  pair  of  observations  at  time  ti  and  tj  in  the  light  curve. 
This  describes  a  DRW  process  with  a  characteristic  damping 
timescale  t  beyond  which  the  structure  function  will  asymptote 
to  a  constant  value  of  SFoo  =  and  a  long-term  standard 
deviation  of  variability  cr  =  V^SF^.  The  short-term  driving 
amplitude  of  variability  a  =  o+JI/t  determines  the  rise  of 
SF(Af)  for  A f  «  r.  We  model  each  individual  light  curve  and 
calculate  a  and  t  ,  along  with  their  likelihood  distributions,  using 
the  method  of  Press  et  al.  (1992),  its  generalization  in  Rybicki 
&  Press  (1992),  and  the  fast  computational  implementation  in 
Rybicki  &  Press  (1995).  The  corresponding  structure  function 
for  the  model  for  time  lag  At  is 

SF(Af)  =  SFoAl  -  e" |A'l/r)1/2.  (3) 

In  fitting  each  individual  light  curve  with  the  DRW  model, 
we  calculate  the  likelihood  of  a  DRW  solution  LDRW,  as  well  as 
the  likelihood  of  a  pure  white  noise  solution  LnoiSe-  Light  curves 
that  are  better  described  by  the  DRW  model  than  white  noise 
at  a  5a  level  will  have  ALnoise  =  ln(LDRW/Lnoise)  >  5.  All  our 
FSRQ  and  BL  Lac  object  light  curves  in  the  Variable  LINEAR 
catalog  fit  the  DRW  model  at  above  this  5a  level.  To  remove 
light  curves  for  which  the  survey  length  is  shorter  than  r  (thus 
leaving  r  unconstrained),  we  also  calculate  the  likelihood  of 
a  runaway  timescale,  L0 0.  Objects  for  which  t  is  constrained 
in  the  DRW  model  will  have  AL^  =  ln(LDRW/Loo)  >  0.05 
(MacLeod  et  al.  2010).  MacLeod  et  al.  (2010)  showed  that  the 


4 


The  Astrophysical  Journal,  760:51  (llpp),  2012  November  20 


Ruan  et  al. 


Figure  4.  Distribution  of  variable  LINEAR  objects  in  SFqo  (asymptotic  value 
of  the  structure  function  on  long  timescales)  and  r  (characteristic  damping 
timescale).  BL  Lac  objects  (red  stars),  FSRQs  (red  circles),  normal  quasars 
(black  crosses),  and  contours  enclosing  90%  (blue)  and  50%  (green)  of  all  other 
objects  are  shown.  The  side  panels  show  the  projected  normalized  histograms 
of  the  distributions  in  SFoo  (top)  and  r  (right)  for  blazars  (BL  Lac  objects  and 
FSRQs,  red  solid),  normal  quasars  (black  dashed),  and  other  objects  (green 
dotted).  Dashed  yellow  lines  of  constant  o  for  o  =  —0.75,  0.0,  0.75  are  also 
shown,  with  o  increasing  toward  the  lower  right. 


DRW  is  an  excellent  fit  to  quasar  light  curves  and  used  cuts 
on  minimum  ALnoiSe  and  minimum  AL^  to  select  quasars  with 
high-quality  light  curves.  In  the  higher  stellar  density  regions 
of  Stripe  82,  MC11  used  cuts  on  maximum  <7,  minimum  r, 
and  minimum  ALnoiSe  to  achieve  impressive  completeness  C  = 
93%  (defined  as  percentage  of  the  total  confirmed  quasars  in  the 
sample  that  are  selected)  and  efficiency  E  >  78%  (percentage  of 
quasar  candidates  selected  that  are  confirmed  quasars)  in  quasar 
selection.  Furthermore,  Kozlowski  et  al.  (2010)  showed  that  the 
method  works  even  in  the  high  stellar  density  regions  of  the 
Magellanic  Cloud  fields. 

We  calculate  the  best-fit  DRW  variability  parameters  of  the 
individual  light  curves  of  the  60  blazars,  155  normal  quasars, 
and  a  random  sample  of  6000  other  objects  (mainly  foreground 
Galactic  stars,  representative  of  the  other  >180,000  objects) 
from  the  Variable  LINEAR  catalog.  To  first  understand  the 
underlying  distribution  of  DRW  variability  parameters  for  these 
different  populations,  we  redshift-correct  our  parameters  as 
outlined  in  Kelly  et  al.  (2009)  using  spectroscopic  redshifts 
from  the  Schneider  et  al.  (2010)  catalog  for  ah  FSRQs  and 
normal  quasars.  For  BL  Lac  objects,  redshifts  are  much  less 
accurate  (and  sometimes  not  possible)  due  to  the  lack  of  strong 
spectral  features.  Since  our  sample  of  BL  Lac  objects  is  drawn 
from  a  variety  of  catalogs,  there  are  sometimes  discordant 
redshifts  reported;  in  such  cases,  we  preferentially  adopt  the 
spectroscopically  derived  redshifts  (or  the  lower  limits)  from 
the  Plotkin  et  al.  (2008)  and  Plotkin  et  al.  (2010)  catalogs.  For 
three  BL  Lac  objects,  there  is  no  redshift  reported  in  any  of 
the  catalogs,  and  so  for  these  three  objects  we  adopt  the  mean 
redshift  found  for  ah  other  LINEAR  BL  Lac  objects,  z  =  0.28. 
While  this  assumption  is  not  ideal,  very  few  extragalactic  objects 
in  the  LINEAR  survey  will  be  at  high  redshifts  due  to  the  shallow 


optical  flux  limit,  so  uncertainties  on  these  three  redshifts  will 
not  strongly  impact  our  results. 

Imposing  the  restrictions  on  ALoo  >0.05  and  ALnoise  >  5 
leaves  us  with  51  blazars  and  121  normal  quasars  in  our  sam¬ 
ple.  Figure  4  shows  the  best-fit  DRW  parameters  r  against  SF^ 
of  blazars,  normal  quasars,  and  other  objects  (mostly  Galac¬ 
tic  stars),  along  with  normalized  histograms  of  their  distri¬ 
butions.  Although  we  use  separate  symbols  for  BL  Lac  ob¬ 
jects  and  FSRQs  in  Figure  4,  we  do  not  separate  these  two 
blazar  sub-populations  in  our  subsequent  analysis  due  to  the 
small  sample  size  and  unknown  biases  in  our  blazar  sam¬ 
ple  stemming  from  possible  correlations  between  variability 
and  physical  properties  of  the  blazars.  Future  surveys  yield¬ 
ing  larger  samples  may  be  able  to  robustly  constrain  differ¬ 
ences  between  these  two  sub-populations  and  provide  further 
insight  into  their  jet  properties.  Figure  4  is  directly  comparable 
to  the  analysis  of  Kozlowski  et  al.  (2010)  and  MA11  and  con¬ 
firms  that  the  structure  functions  of  normal  quasar  light  curves 
tend  to  have  longer  characteristic  timescales  of  variability  r 
and  slightly  larger  SF^  than  stars.  More  importantly,  our  anal¬ 
ysis  also  shows  that  blazars  have  characteristic  timescales  r 
in  between  those  of  normal  quasars  and  other  objects  (mainly 
variable  stars),  as  well  as  larger  values  of  SF^.  The  driving  am¬ 
plitude  on  short  timescales  a  is  related  to  these  two  parameters 
by  the  relation  a  =  SF^/V^7,  as  shown  in  Figure  4.  Figure  5 
shows  the  distribution  for  blazars  and  normal  quasars  in  <7, 
where  normal  quasars  tend  to  peak  at  log  (cr  / mag  yr-1/2)  ~ 
—0.75,  while  blazars  tend  to  peak  at  log  a  ~  0.0.  A  two- 
dimensional  Kolmogorov-Smirnoff  test  between  the  distribu¬ 
tions  of  a  for  blazars  and  normal  quasars  in  Figure  5  gives  a 
/7-value  of  2.04  x  10-5.  For  the  distributions  of  r  in  Figure  4, 
the  /7-value  is  0.032. 

While  MAI  1  performed  quasar  selection  using  this  approach 
by  placing  cuts  on  maximum  a  and  minimum  r ,  we  can  also  see 
from  Figure  4  that  a  minimum  r  cut  would  separate  blazars  from 
other  objects  (mostly  Galactic  stars),  while  a  minimum  a  cut 
can  help  separate  blazars  from  quasars.  However,  it  is  clear  from 
Figure  4  that  highly  efficient  and  complete  quasar  and  blazar 
selection  purely  by  optical  variability  in  the  LINEAR  survey  is 
difficult,  as  there  will  be  either  heavy  contamination  or  a  low 
recovery  rate.  The  main  cause  of  this  is  the  bright  magnitude 
range  (14  <  w linear  <  17)  of  LINEAR,  which  overwhelmingly 
probes  Galactic  stars  and  very  few  actual  AGNs  (e.g.,  there  are 
>180,000  sources  in  the  Variable  LINEAR  catalog,  only  214 
of  which  are  known  normal  quasars  and  blazars),  and  so  deeper 
survey  data  containing  a  substantial  fraction  of  extragalactic 
sources  are  necessary.  However,  currently  available  data  from 
deeper  time-domain  optical  surveys  such  as  SDSS  Stripe  82 
have  comparatively  small  sky  coverage  and  do  not  yield  a  large 
enough  sample  of  the  much  rarer  blazars. 

4.  BLAZAR  SELECTION  IN  FERMI  ERROR  CIRCLES 

Despite  the  overwhelming  stellar  contamination  in  the  bright 
magnitude  range  of  the  LINEAR  survey,  our  proposed  blazar 
selection  method  is  still  useful  with  LINEAR  if  we  restrict 
our  search  to  the  small  regions  of  sky  associated  with  the 
positions  of  Fermi  sources.  This  greatly  reduces  the  stellar 
contamination,  while  simultaneously  providing  a  test  of  a 
novel  method  (based  on  long-term  photometric  variability 
characteristics)  to  associate  ]/-ray  sources  with  counterparts 
at  optical  wavelengths.  We  will  show  the  DRW  parameters 
of  known  Fermi  AGNs  (overwhelmingly  blazars),  quantify 
the  completeness  and  efficiency  of  blazar  selection  using  this 
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Figure  5.  Distributions  of  the  driving  amplitude  on  short  timescales  a  for  blazars 
(solid)  and  normal  quasars  (dashed)  in  the  Variable  LINEAR  catalog. 

method,  and  identify  new  variability- selected  candidate  Fermi 
blazar  counterparts. 

We  first  positionally  match  the  Variable  LINEAR  catalog 
described  in  Section  2.3  to  the  1873  sources  in  the  2FGL  catalog, 
approximating  each  error  ellipse  as  a  circle  with  radius  equal  to 
the  semimajor  axis  of  the  95%  confidence  ellipse.  This  initial 
match  yields  480  variable  LINEAR  objects  in  173  Fermi  error 
circles.  These  480  variable  LINEAR  optical  objects  in  Fermi 
error  circles  are  then  positionally  matched  using  a  3"0  matching 
radius  to  the  929  associated  AGNs  in  the  2nd  Fermi  AGN  catalog 
(Ackermann  et  al.  2011),  restricting  consideration  to  those  that 
have  Bayesian  association  probabilities  >0.8  (this  association 
method  is  described  in  the  Appendix  of  Abdo  et  al.  2010). 
This  yields  a  sample  of  47  known  (i.e.,  confidently  associated) 
Fermi  y -ray-emitting  AGNs.  We  calculate  DRW  parameters  for 
the  LINEAR  light  curves  of  all  480  variable  objects  and  make 
nominals  cut  on  ALnoiSe  >  5  and  ALoo  >  0.05  as  in  Section  3. 

Due  to  the  bright  magnitude  range  of  the  LINEAR  survey,  it 
is  likely  that  the  flux  limit  of  the  LINEAR  survey  is  too  bright 
to  actually  detect  the  optical  counterparts  of  the  faintest  Fermi 
y- ray  sources,  leaving  us  a  population  of  orphan  error  circles 
that  contain  only  LINEAR  contaminants.  Since  these  fainter 
sources  will  also  preferentially  have  larger  error  circles,  we  can 
remove  them  by  considering  further  cuts  on  the  maximum  radii 
of  the  error  circles.  In  Figure  6,  we  show  the  distribution  of 
the  total  number  of  LINEAR  variable  objects  in  Fermi  error 
circles  as  a  function  of  error  circle  radius  r,  for  all  173  error 
circles  that  contain  at  least  one  LINEAR  variable  object,  as 
well  as  the  distribution  for  the  47  error  circles  containing 
a  LINEAR  variable  object  that  matched  to  an  AGN  in  the 
2nd  Fermi  AGN  catalog.  Figure  6  shows  that  the  number  of 
variable  LINEAR  objects  in  error  circles  with  known  y- ray 
AGN  counterparts  peaks  at  2'  and  drops  rapidly;  beyond  10', 
there  is  only  one  Fermi  error  circle  that  contains  a  known  y- 
ray  AGN  counterpart.  However,  the  distribution  of  the  total 
number  of  variable  LINEAR  objects  in  all  error  circles  as  a 
function  of  radius  has  a  long  tail  to  large  r,  as  there  are  a  small 
number  of  large  error  circles  (r  >  15')  with  large  numbers 
of  variable  LINEAR  objects,  as  expected.  This  suggests  that 
a  cut  on  the  error  circle  radius  of  ~10'  is  appropriate,  but  to 
gauge  the  impact  of  such  a  cut,  we  will  consider  in  detail  a 
range  of  8',  10',  and  12'  cuts  on  the  maximum  error  circle  size. 
A  cut  at  r  <  10'  removes  from  consideration  only  one  Fermi 
error  circle  containing  a  known  y- ray  AGN,  while  a  cut  at 
r  <  8'  removes  four.  Scatter  in  the  negative  relation  between 
the  optical  flux  and  y-ray  error  circle  size  of  Fermi  blazars  may 
introduce  some  biases  in  our  resulting  sample  toward  blazars 
with  higher  than  average  optical-to- y-ray  flux  ratios  when  we 
place  a  cut  on  the  error  circle  size.  However,  this  bias  is  expected 
to  be  small,  as  only  ~7%  of  reliably  associated  AGNs  in  the 


Figure  6.  Distribution  of  the  number  of  variable  LINEAR  objects  lying  in  all 
Fermi  error  circles  of  radius  r  (dashed)  and  only  those  lying  in  Fermi  error 
circles  of  radius  r  that  contain  a  variable  LINEAR  object  matched  to  an  AGN 
in  the  2nd  Fermi  AGN  catalog  (solid),  as  a  function  of  r. 

2nd  Fermi  AGN  catalog  with  optical  magnitudes  <17  reside  in 
Fermi  error  circles  >10'  in  radius.  We  note  that  although  the 
Fermi  AGN  association  success  rate  of  Ackermann  et  al.  (201 1) 
as  a  function  of  r  will  also  affect  Figure  6,  this  effect  is  secondary 
to  the  LINEAR  flux  limit  for  the  high  association  probabilities 
considered  here. 

Figure  7  (left  panel)  shows  the  distribution  in  r  and  SFoo  of 
LINEAR  optical  counterparts  to  Fermi  y- ray  AGNs,  as  well 
as  all  other  variable  objects  in  Fermi  2FGL  error  circles  of 
r  <  10'.  We  do  not  apply  cosmological  redshift  corrections 
here  for  blazar  selection.  Similar  to  the  conclusions  drawn 
from  Figure  4,  there  is  clear  separation  between  y  -ray-emitting 
AGNs  and  other  LINEAR  variables  in  DRW  parameter  space. 
The  significant  normal  quasar  population  seen  in  Figure  4 
is  minimized  when  looking  only  in  Fermi  error  circles.  As 
expected,  almost  all  y-ray  AGNs  recovered  in  our  sample  have 
blazar-like  variability.  Figure  7  (right  panel)  shows  that  variable 
objects  in  error  circles  with  radius  r  >  10'  largely  have  SF^ 
and  r  similar  to  that  of  stars  and  are  thus  mostly  contaminants, 
rather  than  new  blazar  candidates  that  were  removed  by  the  cut 
on  r  <  10'. 

We  can  select  y- ray  AGNs  by  placing  cuts  on  minimum  r 
and  a.  For  different  log  r  and  log  <7,  we  estimate  efficiency  and 
completeness  of  y- ray  AGN  selection.  The  results  for  a  10'  cut 
on  the  error  circles  are  shown  in  Figure  8.  Selection  cuts  at  large 
minimum  t  and  minimum  a  yield  efficiencies  of  0  as  there  are 
no  such  objects.  Figure  9  shows  the  maximum  completeness 
achievable  by  tweaking  the  selection  cuts  on  minimum  r  and 
cr,  as  a  function  of  efficiency.  This  shows  how  the  completeness 
changes  as  a  function  of  efficiency  in  y- ray  AGN  selection, 
for  error  circles  with  a  maximum  radius  of  8',  10',  and  122  As 
expected,  the  efficiency  decreases  as  the  maximum  size  of  error 
circles  considered  increases,  at  a  set  completeness. 

By  jointly  maximizing  efficiency  and  completeness  via  the 
quantity  VE2  +  C2  from  the  curves  in  Figure  9,  we  can  achieve 
E  ^  88%  and  C  =  88%  for  Fermi  error  circles  of  r  <  8'  using 
cuts  on  log  r  >  1.35  and  log  a  >  —1.20.  For  error  circles  of 
r  <  10',  E  ^  76%  and  C  =  86%  using  cuts  on  log  r  >  1.53  and 
log  <7  >  —0.90,  and  for  error  circles  of  r  <  12',  E  ^  70%  and 
C  =  86%  using  cuts  on  log  t  >  1.53  and  log  o  >  —0.90.  This 
demonstration  verifies  that  we  are  able  to  select  (in  this  case, 
recover  in  double-blind  fashion)  y  -ray-emitting  AGNs  in  Fermi 
error  circles  with  high  completeness  and  efficiency  by  modeling 
their  optical  light  curves  as  a  DRW  and  imposing  selection  cuts 
on  the  variability  parameters. 

Our  definition  for  completeness  is  similar  to  the  convention 
of  many  other  recent  time-domain  AGN  studies  (e.g.,  MA11). 
The  completeness  we  calculate  here  is  specifically  for  the  se¬ 
lection  of  variable  LINEAR  objects  matched  to  Fermi  AGNs  in 
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Figure  7.  Left:  distribution  of  variable  LINEAR  objects  in  Fermi  error  circles  of  radius  r  <  10',  in  SF^  and  r,  for  reliably  associated  Fermi  AGNs  (red  stars)  and 
other  objects  (black  crosses).  The  underlying  contours  from  Figure  4  enclosing  90%  (blue)  and  50%  (green)  of  all  other  objects  in  the  Variable  LINEAR  catalog  are 
shown  in  the  background  for  reference.  Right:  same  as  the  left  panel,  but  for  Fermi  error  circles  >  10'  in  radius. 


min  log((j/(mag/yr1/2)) 


min  log(<T/(mag/yr1/2)) 


Figure  8.  Completeness  and  efficiency  of  selection  of  confirmed  Fermi  AGNs  in  error  circles  of  radius  r  <  10',  as  a  function  of  selection  cuts  on  minimum  r  and  a. 


Ackermann  et  al.  (2011)  with  Bayesian  association  probabili¬ 
ties  >0.8,  and  the  efficiency  is  calculated  assuming  that  these 
associations  are  all  correct.  However,  this  efficiency  is  a  lower 
bound,  as  it  will  increase  if  some  of  the  contaminants  we  en¬ 
counter  in  selection  are  actually  blazars.  This  may  occur  in  our 
analysis  if  some  LINEAR  variables  that  pass  all  DRW  blazar 
selection  cuts  and  did  not  match  to  an  AGN  in  the  2nd  Fermi 
AGN  catalog  may  actually  be  additional  y- ray  AGN  counter¬ 
parts,  not  already  recognized  as  such.  Indeed,  this  is  likely  the 
case  for  at  least  one  object,  lying  in  the  error  circle  of  Fermi 
source  2FGL  J1649.6+5238. 

In  the  2nd  Fermi  AGN  catalog,  Fermi  source  2FGL 
J 1649.6+5238  is  listed  as  associated  with  the  radio  source 
87GB  164812.2+524023  at  a  0.0  Bayesian  probability,  calcu¬ 
lated  based  on  the  local  density  of  sources  from  catalogs  of 
likely  counterparts.  This  radio  source  is  positionally  coinci¬ 
dent  with  a  LINEAR  object  at  R.A.  =  16h49m24s.99,  decl.  = 
52°35/15"05,  which  has  DRW  parameters  calculated  from  its 
LINEAR  light  curve  of  log  r  =  1 .936,  log  a  =  —0.026,  and  log 
SFqo  =  0.942,  very  similar  to  that  of  confirmed  y- ray  AGNs  in 
Figure  7.  This  LINEAR  source  was  not  counted  as  a  y-ray 


AGN  in  our  calculations  of  efficiency  and  completeness  be¬ 
cause  its  Bayesian  probability  of  association  was  below  0.8  in 
Ackermann  et  al.  (2011).  Furthermore,  this  LINEAR  counter¬ 
part  also  did  not  match  to  any  known  blazar  in  the  BL  Lac  object 
and  quasar  catalogs  used  in  Section  2.2. 

Based  on  these  blazar-like  DRW  parameters  calculated  from 
the  light  curve,  we  suggest  this  LINEAR  variable  object 
as  the  plausible  optical  counterpart  to  Fermi  source  2FGL 
J 1649. 6+5 23,  as  well  as  increased  confidence  in  associating 
87GB  164812.2+524023  as  the  radio  counterpart.  We  note  that 
this  Fermi  source  is  also  associated  with  this  radio  source  with 
a  0.82  probability  using  the  log  A/-log  S  method  in  Ackermann 
et  al.  (2011),  based  on  observed  properties  of  candidate  radio 
counterparts.  The  origin  of  this  disparity  lies  in  the  differences  in 
the  approaches  used  to  calculate  the  Bayesian  and  log  A-log  S 
(radio)  association  probabilities.  We  opted  to  use  the  Bayesian 
probability  >0.8  criterion  for  known  y-ray  AGNs,  but  the  use 
of  the  log  A/-log  S  probability  also  recovers  this  radio  source 
as  a  likely  y- ray  AGN.  It  is  not  our  intention  to  compare  the 
two  statistical  approaches,  but  rather  to  point  out  that  the  dis¬ 
crepancies  in  the  results  of  these  two  methods  may  lead  to 
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Table  1 

J2000  Positions  of  New  Candidate  LINEAR  Blazars  in  Fermi  Error  Circles  with  r  <  10'  Selected  by  Variability 


R.A. 

(deg) 

Decl. 

(deg) 

Fermi  2FGL  Source 

log 

(r/days) 

log 

(<r/(mag  yr-1/2)) 

2nd  Fermi  AGN  Catalog  Association 

128.5027 

42.2570 

2FGL  J0834.3+4221 

1.89 

-0.20 

OJ  45 

153.3687 

34.4883 

2FGL  J1 01 3.6+3434 

1.58 

-0.28 

Unassociated 

172.2141 

-5.6568 

2FGLJ1 129.0-0532 

2.20 

0.02 

Unassociated 

189.5802 

-20.0237 

2FGLJ1238. 1-1953 

2.07 

-0.11 

PMN  J1238— 1959 

219.1272 

23.2824 

2FGL  J1436.9+2319 

1.66 

-0.29 

PKS  B 1434+235 

220.4843 

43.9407 

2FGL  J1442.0+4352 

2.00 

-0.26 

BZB  J1442+4348 

232.8992 

57.4483 

2FGL  J153 1.0+5725 

2.50 

-0.04 

Unassociated 

232.9009 

57.3541 

2FGL  J153 1.0+5725 

2.36 

-0.58 

Unassociated 

246.9883 

32.4363 

2FGL  J1627. 8+3219 

2.30 

-0.34 

Unassociated 

250.2259 

11.7870 

2FGL  J164 1.0+1 141 

1.89 

-0.20 

MG1  J164058+1144 

252.3542 

52.5875 

2FGL  J1649.6+5238 

1.93 

-0.02 

Unassociated 

331.7412 

-0.4681 

2FGL  J2206.6— 0029 

1.80 

0.27 

PMN  J2206— 0031 

efficiency 

Figure  9.  Maximum  completeness  achievable  for  a  given  efficiency  of  selection, 
for  confirmed  Fermi  AGNs  in  error  circles  of  radius  r  <  8'  (dashed),  r  <  10' 
(solid),  and  r  <  12'  (dotted)  as  a  function  of  efficiency. 

different  results  in  analysis.  In  cases  where  these  two  methods 
lead  to  highly  discrepant  results,  our  completely  independent  op¬ 
tical  variability-based  method  can  provide  valuable  additional 
information. 

Finally,  we  note  that  our  variability-based  approach  may 
be  further  applicable  to  mass  identification  of  unassociated 
y-ray  sources  using  deeper  time-domain  optical  surveys.  Data 
from  current  surveys  such  as  the  PTF,  with  single-epoch 
depth  of  r  ~  20.6  mag,  30  epochs  of  observation  over 
8200  deg2  of  sky,  and  ~0.01  mag  repeatability,  can  be  used 
to  identify  unassociated  Fermi  blazars  using  the  method  we 
have  presented  here.  Future  surveys  such  as  LSST,  with  single¬ 
epoch  depth  of  r  ~  24.5  mag  over  20,000  deg2  of  sky  and  ~103 
epochs  per  source,  can  be  used  to  vastly  increase  the  existing 
sample  size. 

5.  DISCUSSION 
5.1.  Efficiency  and  Completeness 

In  Section  3,  we  noted  that  blazar  selection  using  the  Variable 
LINEAR  catalog  across  the  full  field  of  the  survey  (as  contrasted 
to  blazar  selection  only  within  Fermi  error  circles)  was  not 
efficient  and  complete.  At  first,  this  may  seem  to  be  in  stark 
contrast  to  the  success  of  quasar  variability  selection  using 


a  similar  DRW  modeling  approach  (e.g.,  MA11;  Butler  & 
Bloom  2011;  Kozlowski  et  al.  2010).  This  is  primarily  a 
consequence  of  the  shallow  depth  of  the  LINEAR  survey  rather 
than  a  reflection  of  the  ultimate  performance  expected  from 
our  approach  (especially  for  fainter  blazars).  At  the  depth  of 
SDSS  Stripe  82  (flux  limit  of  g  ~  20.5),  objects  with  low- 
redshift  quasar-like  colors  account  for  ~63%  of  all  variable 
objects,  with  Galactic  stars  making  up  the  vast  majority  of  the 
remaining  variable  objects  (Sesar  et  al.  2007).  In  contrast,  at  the 
shallow  depth  of  the  LINEAR  survey  (flux  limit  of  r  ~  17), 
only  ~0.08%  of  objects  in  the  Variable  LINEAR  catalog  are 
quasars  (see  Section  2.3).  As  such,  quasar  selection  by  imposing 
appropriate  maximum  d  and  minimum  x  cuts  (as  done  in  MAI  1) 
on  the  LINEAR  quasar  sample  in  Figure  4  will  not  yield  E 
beyond  a  few  percent  for  any  value  of  C  >  50%  due  to  scatter  in 
the  parameters  of  the  enormously  larger  population  of  variable 
stars.  Blazar  selection  in  LINEAR  will  be  similarly  inefficient. 
These  issues  are  exacerbated  by  the  large  errors  on  the  estimated 
DRW  parameters  as  compared  to  SDSS,  caused  by  the  larger 
photometric  uncertainties  in  LINEAR.  In  any  case,  current 
and  future  time-domain  surveys  such  as  PTF  and  LSST  will 
alleviate  both  these  issues  by  going  many  magnitudes  deeper 
while  providing  ~1%  photometry. 

As  noted  in  Section  4,  our  calculated  efficiencies  of  blazar 
selection  are  lower  limits,  since  calculating  the  true  effi¬ 
ciency  requires  correct  identification  of  every  object  as  either  a 
blazar  or  contaminant.  There  are  many  contaminating  variable 
LINEAR  objects  that  fall  in  the  Fermi  error  circles  we  con¬ 
sidered  in  our  analysis  in  Section  4  that  have  blazar-like  vari¬ 
ability  (i.e.,  have  DRW  parameters  similar  to  blazars  and  are 
counted  as  contaminants  when  jointly  optimizing  complete¬ 
ness  and  efficiency  for  known  y- ray  blazar  selection).  Some 
of  these  contaminants  are  certainly  blazars  and  may  actually 
be  the  y- ray  source.  Securing  identification  of  these  blazar-like 
variable  objects  will  require  spectroscopic  follow-up.  We  list 
in  Table  1  the  positions,  parent  2FGL  error  circle,  and  DRW 
parameters  of  12  such  new  candidate  LINEAR  variable  ob¬ 
jects  with  blazar-like  optical  variability  that  fall  in  Fermi  2FGL 
error  circles  with  r  <  10'.  We  also  list  any  previously  associ¬ 
ated  counterparts  of  each  parent  2FGL  source  from  Ackermann 
et  al.  (2011)  with  Bayesian  association  probabilities  >0.8.  The 
LINEAR  blazar  candidates  are  not  positionally  coincident  with 
these  previously  associated  counterparts;  rather,  they  are  addi¬ 
tional  y-ray  blazar  candidates  in  their  respective  error  circles. 
Since  it  is  possible  that  some  Fermi  error  circles  may  contain 
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more  than  one  y  -ray-emitting  object,  these  LINEAR  blazar  can¬ 
didates  are  still  worthy  of  spectroscopic  follow-up. 

We  further  note  that  the  efficiencies  we  have  calculated  in 
y- ray  blazar  selection  in  Section  4  are  conservative,  as  we  have 
considered  all  variable  LINEAR  objects  lying  in  all  Fermi  error 
circles  below  a  certain  size.  However,  ~10%  of  sources  in  the 
2FGL  are  associated  with  non-AGN  objects  (pulsars,  supernova 
remnants,  etc.),  while  ~30%  of  sources  are  unassociated.  We 
chose  to  include  all  these  error  circles  in  our  analysis  (and  not 
just  those  containing  known  Fermi  blazars)  in  a  double-blind 
test,  as  it  is  a  more  faithful  demonstration  of  the  efficacy  of  this 
approach. 

We  would  like  to  compare  the  efficiency  and  completeness 
of  our  variability-based  Fermi  blazar  selection  method  to  other 
methods.  However,  a  direct  comparison  is  difficult;  few  stud¬ 
ies  have  attempted  to  use  characteristic  properties  of  blazars 
to  systematically  select  Fermi  blazars  in  double-blind  fash¬ 
ion  as  we  have  done.  Furthermore,  the  completeness  and  ef¬ 
ficiency  of  those  methods  are  almost  never  jointly  optimized 
(as  done  in  Section  4).  For  example,  Massaro  et  al.  (2011)  used 
mid-IR  color  selection  to  recover  known  blazars  in  the  WISE 
survey  and  extended  their  approach  to  find  new  Fermi  blazar  can¬ 
didates  (Massaro  et  al.  2012).  By  parameterizing  the  similarity 
of  the  mid-IR  colors  of  WISE  sources  to  W/SF-detected  Fermi 
blazars,  Massaro  et  al.  (2012)  imposed  a  mid-IR  color-based  se¬ 
lection  criterion  that  is  able  to  select  known  Fermi  blazars  with 
87%  completeness  over  the  full  region  of  sky  surveyed  by  WISE 
in  its  first  year.  However,  Massaro  et  al.  (2012)  did  not  jointly 
calculate  the  completeness  and  efficiency  of  their  approach. 
Their  mid-IR  color  selection  method  has  the  advantage  of  the 
all-sky  coverage  of  the  WISE  survey.  However,  mid-IR  color 
selection  is  dependent  on  the  degree  to  which  non-thermal  jet 
emission  dominates  the  mid-IR  emission  of  individual  blazars. 
Indeed,  Massaro  et  al.  (2011)  find  that  FSRQs  generally  have 
WISE  colors  closer  to  normal  quasars  than  to  BL  Lac  objects. 
This  may  be  due  to  significant  thermal  mid-IR  emission  in  the 
SEDs  of  FSRQs,  similar  to  that  in  normal  quasars  (e.g.,  from 
dust  emission;  see  Plotkin  et  al.  2012),  thus  making  blazar  se¬ 
lection  more  difficult. 

Massaro  et  al.  (2012)  apply  their  selection  method  to  WISE 
sources  in  unidentified  Fermi  error  circles  to  find  297  y- ray 
blazar  candidates  in  156  out  of  313  unidentified  Fermi  error  cir¬ 
cles.  However,  their  estimates  of  the  efficiency  of  their  method 
by  systematically  offsetting  the  position  of  each  Fermi  error  cir¬ 
cle  to  assess  the  number  of  random  associations  gave  262  false 
blazar  candidates  (using  the  same  mid-IR  color  selection  crite¬ 
ria),  suggesting  that  the  efficiency  of  this  approach  may  be  low. 
Nevertheless,  mid-IR  color-  and  variability-based  selections  of 
blazars  are  independent  and  complementary  methods,  with  dif¬ 
ferent  selection  biases.  A  combination  of  both  approaches,  using 
the  full  WISE  survey  as  well  as  data  from  current  time-domain 
optical  imaging  surveys,  may  be  highly  efficient  and  complete. 

5.2.  Implications  for  Blazar  Variability 

Our  finding  that  optical  blazar  variability  is  well  described  as 
a  DRW  process  with  distinct  variability  parameters  may  have 
interesting  potential  implications  for  blazar  jet  physics.  The  SFoo 
and  t  DRW  parameters  should  both  be  affected  by  relativistic 
effects.  For  example,  the  observed  blazar  characteristic  damping 
timescale  Tbiz,0bs  (after  correcting  for  cosmological  redshift) 
should  be  shortened  in  comparison  to  the  interval  Tbiz,rest  in  the 
rest  frame  of  the  emitting  material  due  to  the  relativistic  motion, 
with  Tbiz,rest  =  ^b\z,obs,  where  8  is  the  kinematic  Doppler  factor. 
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Figure  10.  Top:  distribution  of  calculated  kinematic  Doppler  factors  of  blazars 
in  the  Variable  LINEAR  catalog,  assuming  that  the  rest-frame  characteristic 
timescale  of  variability  for  each  blazar  is  equal  to  that  of  LINEAR  normal 
quasars  (320  days).  Bottom:  distribution  of  calculated  rest-frame  characteristic 
timescale  of  variability  for  blazars  in  the  Variable  LINEAR  catalog,  assuming 
a  kinematic  Doppler  factor  <5  =  10. 


While  Tbiz,rest  is  uncertain  and  difficult  to  measure,  we  can  place 
constraints  on  it  through  comparisons  to  rqso,  the  characteristic 
DRW  timescale  of  normal  quasars  (which  are  not  affected  by 
relativistic  effects),  and  independent  measurements  of  8. 

In  Figure  4,  the  distribution  of  r  peaks  at  Tbiz,0bs  ~  80  days 
for  blazars  and  rqso  ~  320  days  for  normal  quasars.  While 
these  distributions  have  large  tails,  their  relative  peaks  pro¬ 
vide  an  adequate  basis  for  comparison.  If  the  underlying  vari¬ 
ability  in  quasars  ultimately  also  drives  blazar  variability,  then 
8  ~  Tqso/ibiz,obs  ~  4.  Figure  10  shows  the  distribution  of  this 
estimated  8  for  blazars  in  the  Variable  LINEAR  catalog,  assum¬ 
ing  rqso  =  320  days.  This  is  smaller  than  estimates  of  the  typical 
blazar  Doppler  factor  from  radio  observations  of  Fermi- detected 
blazars  (5  ~  10-30;  Savolainen  et  al.  2010)  and  may  imply  that 
the  variability  of  the  non-thermal  jet  emission  in  the  jet  rest 
frame  (e.g.,  due  to  shocks  in  the  jet)  occurs  on  longer  timescales 
than  the  thermal  disk  variability  of  normal  quasars. 

We  can  instead  adopt  the  observed  typical  Doppler  factors  for 
blazars  and  calculate  Tbiz,rest  for  each  individual  blazar.  Ideally, 
we  would  use  a  measured  value  of  8  for  each  individual  blazar  in 
this  calculation,  but  direct  measurements  of  8  from  the  literature 
are  sparse  and  inhomogeneous.  In  Figure  10  (bottom  panel),  we 
show  the  distribution  of  rarest  assuming  a  reasonable  8  =  10 
for  all  blazars.  The  peak  in  Figure  10  occurs  at  rarest  ~  3  years, 
which  is  longer  than  the  rqso  ~  320  days  we  measure  for 
LINEAR  normal  quasars.  The  discrepancy  between  8  for  Fermi 
blazars  implied  by  radio  observations  and  8  calculated  assuming 
intrinsic  variability  similar  to  normal  quasars  suggests  different 
stochastic  mechanisms  driving  the  variability  in  the  disk  and  the 
jet.  While  a  detailed  physical  interpretation  of  this  is  currently 
unclear,  our  long-term  optical  variability  results  may  provide 
additional  constraints  on  models  of  blazar  jets. 

5.3.  Short-timescale  Variability  of  Quasars  and  Blazars 

Observations  of  intraday  variability  (also  referred  to  as 
“microvariability”)  of  normal  quasars  have  shown  a  puzzling 
diversity  of  properties.  Not  all  quasars  exhibit  this  phenomenon, 
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but  some  show  variability  on  the  order  of  ~0.01  mag  on  ~1  day 
timescales  (Gopal-Krishna  et  al.  2003;  Stalin  et  al.  2004,  2005; 
Gupta  &  Joshi  2005;  Mushotzky  et  al.  201 1).  The  cause  of  quasar 
variability  on  such  short  timescales  is  an  unresolved  problem, 
although  the  presence  of  strong  intraday  variability  is  correlated 
with  radio-loudness,  and  it  has  long  been  suspected  that  a 
weak  jet  component  may  be  both  the  source  of  radio  emission 
and  the  cause  of  the  stronger  variability.  Kelly  et  al.  (2009) 
showed  that  the  DRW  model  successfully  predicts  amplitudes 
of  variability  of  <0.02  mag  over  ~8  hr  for  quasars,  consistent 
with  observations  from  quasar  monitoring.  This  successful 
prediction  of  the  short-timescale  variability  is  noteworthy  since 
the  sampling  intervals  of  the  quasar  light  curves  used  in  Kelly 
et  al.  (2009)  are  all  >2  days  and  thus  do  not  actually  sample 
such  intraday  timescales. 

From  the  log  a  distributions  for  normal  quasars  and  blazars  in 
Figure  5,  we  can  calculate  the  standard  deviation  of  the  expected 
variability  on  ~1  day  timescales,  approximated  as  &VA t.  For 
quasars,  the  peak  of  the  distribution  occurs  at  ~0.2  mag  yr-1/2, 
thus  predicting  intraday  variability  ~0.01  mag  on  1  day 
timescales.  For  blazars,  the  peak  occurs  at  ~1  mag  yr_1/2,  thus 
predicting  variability  of  ~0.05  mag  on  1  day  timescales  in  the 
observed  frame.  Previous  studies  of  blazar  micro  variability  have 
often  focused  on  the  most  variable  blazars  and  have  reported 
night- to-night  variations  as  high  as  1.0  mag  (e.g.,  Carini  et  al. 
1991;  Ghosh  et  al.  2000).  The  intraday  variability  we  have  cal¬ 
culated  for  LINEAR  blazars  is  from  an  untargeted  (i.e.,  less 
biased)  survey  and  is  likely  more  representative  of  the  intraday 
variability  of  blazars  as  a  whole. 

Despite  the  large  photometric  uncertainties  in  the  LINEAR 
survey,  there  are  a  significant  number  of  normal  quasars  and 
blazars  in  Figure  5  that  do  have  extreme  levels  of  short- 
timescale  variability,  above  the  photometric  uncertainty.  Nearly 
all  LINEAR  objects  have  photometric  uncertainties  below 
0.1  mag  (and  reaching  as  low  as  ^0.04)  at  <17  mag  (see 
Figure  2).  Conservatively  assuming  that  variability  on  1  day 
timescales  needs  to  be  >0.1  mag  in  order  to  be  detectable 
above  the  photometric  uncertainty  implies  a  log  a  >  0.28.  In 
Figure  5,  18  of  119  normal  quasars  and  10  of  51  blazars  in 
our  sample  are  variable  above  this  criterion  (and  both  these 
fractions  would  increase  if  we  more  carefully  considered  the 
photometric  uncertainty  limit  of  each  individual  object).  Thus, 
for  these  objects,  the  intraday  variability  is  above  the  LINEAR 
photometric  uncertainty,  and  the  >0.1  mag  intraday  variability 
is  significantly  above  most  previous  observations  of  normal 
quasars.  The  source  of  this  extreme  variability  for  a  subset  of 
quasars  is  unknown  and  warrants  additional  investigation. 

6.  CONCLUSIONS 

We  have  used  light  curves  from  the  LINEAR  optical  imaging 
survey  to  study  the  variability  of  blazars  and  to  compare  their 
variability  to  that  of  normal  quasars.  In  contrast  to  many 
earlier  studies  of  blazar  variability,  this  study  is  the  first  to 
use  a  reasonably  large,  homogeneous  sample  of  individual 
light  curves  from  a  wide-field  survey,  where  the  distribution 
of  various  variability  properties  for  blazars  and  quasars  can  be 
calculated  and  compared. 

In  a  time-domain  survey,  blazars  are  among  the  most  violently 
variable  objects  detected  and  are  generally  more  variable  than 
normal  quasars.  We  estimate  that  ~36%  of  blazars  but  only 
~  1 1  %  of  quasars  are  variable  at  the  extreme  >0.15  mag  intrinsic 
rms  variability  level.  We  show  that,  like  quasars,  blazar  light 


curves  are  well  fit  by  the  DRW  model  for  variability,  but  the 
blazars  lie  in  distinct  variability  parameter  space  with  higher  t 
than  stars  and  higher  a  than  normal  quasars.  This  suggests  that 
blazars  can  be  selected  with  high  efficiency  E  and  completeness 
C  by  imposing  minimum  selection  cuts  on  r  and  a . 

Due  to  the  overwhelming  numbers  of  variable  stars  in  the 
bright  magnitude  regime  probed  by  LINEAR,  it  is  difficult  to 
select  blazars  (and  quasars)  at  high  efficiency  and  completeness 
in  the  full  LINEAR  sample  using  optical  variability  methods 
alone.  We  instead  examine  a  more  focused  test  to  select  y  -ray- 
emitting  blazars  in  Fermi  error  ellipses  from  the  2FGL  catalog. 
We  calculate  DRW  parameters  for  light  curves  of  all  LINEAR 
variable  objects  within  Fermi  error  ellipses  and  place  minimum 
selection  cuts  on  r  and  a.  Using  this  approach,  we  are  able  to 
recover  the  corresponding  y  -ray-emitting  AGN  counterparts  in 
the  2nd  Fermi  AGN  catalog  with  E  ^  88%  and  C  =  88%  for 
95%  confidence  error  ellipses  with  semimajor  axis  r  <  8'  and 
E  ^  70%  and  C  =  86%  for  r  <  12'. 

Our  y- ray  blazar  selection  has  uncovered  a  variable  LINEAR 
object,  coincident  with  radio  source  87GB  164812.2+524023, 
in  the  error  ellipse  of  Fermi  source  2FGL  J1649. 6+5238.  This 
object  was  not  associated  with  the  Fermi  source  in  the  2nd  Fermi 
AGN  catalog  using  Bayesian  probabilities.  Our  analysis  shows 
that  this  object  has  optical  variability  properties  consistent  with 
y  -ray-emitting  blazars  and  is  likely  to  be  the  y- ray  source.  We 
find  a  total  of  12  objects  with  LINEAR  variability  parameters 
similar  to  blazars  lying  in  Fermi  error  circles  with  r  <  10' 
(see  Table  1).  Confirmation  of  these  variability-selected  blazar 
candidates  will  require  spectroscopic  follow-up. 

Our  results  suggest  that  the  variability  of  the  non-thermal 
jet  emission  in  blazars  is  stochastic  in  nature,  with  higher 
amplitudes  and  shorter  observed  timescales  in  comparison  to 
the  variability  of  normal  quasars,  likely  due  to  the  effects  of 
relativistic  beaming.  Assuming  a  reasonable  Doppler  factor  of 
10,  we  estimate  that  blazars  are  characteristically  variable  on 
timescales  of  ~3  years  in  the  rest  frame  of  the  emitting  material, 
longer  than  the  ~320  day  characteristic  disk  flux  timescale 
for  quasars.  This  suggests  that  different  physical  mechanisms 
dominate  the  observed  variability  from  blazars  and  quasars, 
likely  connected  to  the  jet  and  to  the  disk,  respectively.  The 
fitted  parameters  imply  that  blazars  have  a  typical  intraday 
variability  amplitude  of  ~0.05  mag,  compared  to  ~0.01  mag 
for  normal  quasars.  Furthermore,  there  is  a  significant  fraction 
(~15%)  of  normal  quasars  that  exhibit  large  intraday  variability 
of  >0.1  mag,  detectable  above  the  photometric  uncertainty. 
The  source  of  the  extreme  variability  of  these  quasars  is 
unclear. 

We  argue  that  variability-based  blazar  selection  is  likely  to 
be  highly  efficient  and  complete  in  deeper  optical  time-domain 
imaging  surveys,  and  that  the  variability-based  blazar  selection 
method  presented  in  this  paper  is  capable  of  greatly  increasing 
the  number  of  known  blazars.  Our  approach,  based  on  long¬ 
term  photometric  variability  characteristics  of  a  large  sample  of 
individual  blazar  light  curves,  may  also  bring  a  new  perspective 
on  accretion  disk  and  jet  physics.  Finally,  this  work  typifies 
the  fruitful  ancillary  science  made  possible  by  combining  data 
from  surveys  as  diverse  in  all  respects  as  LINEAR,  SDSS,  and 
Fermi. 
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